﻿/*
    Copyright 2010, 2011 Eric Wong
	contact@optimalcycling.com
	http://www.optimalcycling.com
	Optimal Cycling - Advanced power pacing program for cyclists
	
	This file is part of Optimal Cycling.

    Optimal Cycling is free software: you can redistribute it and/or modify
    it under the terms of the GNU Lesser General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Optimal Cycling is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU Lesser General Public License for more details.

    You should have received a copy of the GNU Lesser General Public License
    along with Optimal Cycling.  If not, see <http://www.gnu.org/licenses/>.
*/

#ifndef _FILE_INTERPOLATE_H_
#define _FILE_INTERPOLATE_H_

#ifndef __OPENCL_VERSION__
	#include "Options.h"
	#include "Constants.h"
	#include "Vector.h"
#endif

// Note: Do not use native_divide for the interpolation functions. It is much slower than dividing normally.

/// <summary>
/// Efficiently searches throw an array of independent values to interpolate the dependent values.
/// </summary>
/// <param name="valFind">The value we want to find in the independentData.</param>
/// <param name="independentData">The independent values. Must have at least 1 value and be sorted in ascending order.
///		Must be equal in size to dependentData.</param>
/// <param name="dependentData">The dependent values. Must have at least 1 value and be sorted in ascending order.
///		Must be equal in size to independentData.</param>
/// <param name="dataCount">The number of values in each of independentData and dependentData.</param>
/// <returns>Returns dependentData[0] if valFind is below all independentData.
///		Returns dependentData[dataCount - 1] if valFind is above all independentData.
///		Otherwise, returns the interpolated value of dependentData determined by finding valFind in independentData.</returns>
inline real InterpolateBinary(const real valFind, __global const real* independentData, __global const real* dependentData, 
							  const sizeT dataCount)
{
    sizeT left = 0;
    sizeT right = dataCount - 1;
	real valDiff;
	sizeT indexer;

    if (independentData[0] >= valFind)
    {
        return dependentData[0];
    }
    else if (independentData[dataCount - 1] <= valFind)
    {
        return dependentData[dataCount - 1];
    }

    while (right - left > 1)
    {
        indexer = left + ((right - left) / 2);

        if (independentData[indexer] < valFind)
        {
            left = indexer;
        }
        else
        {
            right = indexer;
        }
    }

	valDiff = independentData[right] - independentData[left];
    if (valDiff < SMALL_REAL)
        return dependentData[left];
    else
        return mix(dependentData[left], dependentData[right], (valFind - independentData[left]) / valDiff);
}


// In real2, the first val is the independent and the second val is dependent
inline real InterpolateBinaryReal2(const real valFind, __global const real2* data, const sizeT dataCount)
{
    sizeT left = 0;
    sizeT right = dataCount - 1;
	real valDiff;
	sizeT indexer;

    if (data[0].s0 >= valFind)
    {
        return data[0].s1;
    }
    else if (data[dataCount - 1].s0 <= valFind)
    {
        return data[dataCount - 1].s1;
    }

    while (right - left > 1)
    {
        indexer = left + ((right - left) / 2);

        if (data[indexer].s0 < valFind)
        {
            left = indexer;
        }
        else
        {
            right = indexer;
        }
    }

	valDiff = data[right].s0 - data[left].s0;
    if (valDiff < SMALL_REAL)
        return data[left].s1;
    else
        return mix(data[left].s1, data[right].s1, (valFind - data[left].s0) / valDiff);
}

///// <summary>
///// Efficiently searches throw an array of independent values to interpolate the dependent values.
///// </summary>
///// <param name="valFind">The value we want to find in the independentData.</param>
///// <param name="independentData">The independent values residing in a linear array with independent values of other data sets. 
/////		Must have at least 1 value and be sorted in ascending order.
/////		Number of independent data values we can use must be equal in size to dependentData.</param>
///// <param name="startIndex">The index at which to start look for data points in independentData.</param>
///// <param name="dependentData">The dependent values. Must have at least 1 value and be sorted in ascending order.</param>
///// <param name="dataCount">The number of values in each of independentData the number of values we can use in dependentData.</param>
///// <returns>Returns dependentData[startIndex] if valFind is below all independentData.
/////		Returns dependentData[startIndex + dataCount - 1] if valFind is above all independentData.
/////		Otherwise, returns the interpolated value of dependentData determined by finding valFind in independentData.</returns>
//inline real InterpolateBinaryExReal(const real valFind, __global const real * independentData, const sizeT startIndex,
//							  __global const real * dependentData, const sizeT dataCount);



inline void InterpolateBinaryParams(const real valFind, __global const real* independentData, const sizeT dataCount,
									 sizeT* left, sizeT* right, real* valDiff, real* controlVal)
{
    if (independentData[0] >= valFind)
    {
		*left = 0;
		*right = 0;
		*valDiff = 0;
		*controlVal = 0;
        return;
    }
    else if (independentData[dataCount - 1] <= valFind)
    {
		*left = dataCount - 1;
		*right = dataCount - 1;
		*valDiff = 0;
		*controlVal = 0;
        return;
    }

    *left = 0;
    *right = dataCount - 1;

    while (*right - *left > 1)
    {
        sizeT indexer = *left + ((*right - *left) / 2);

        if (independentData[indexer] < valFind)
            *left = indexer;
        else
            *right = indexer;
    }

	*valDiff = independentData[*right] - independentData[*left];
    if (!(*valDiff < SMALL_REAL))
		*controlVal = (valFind - independentData[*left]) / *valDiff;
	else
		*controlVal = (real)0;
}

#endif		// _FILE_INTERPOLATE_H_


